Exploring benzoxaborole derivatives as carbonic anhydrase inhibitors: a structural and computational analysis reveals their conformational variability as a tool to increase enzyme selectivity

Abstract Recent studies identified the benzoxaborole moiety as a new zinc-binding group able to interact with carbonic anhydrase (CA) active site. Here, we report a structural analysis of benzoxaboroles containing urea/thiourea groups, showing that these molecules are very versatile since they can bind the enzyme assuming different binding conformations and coordination geometries of the catalytic zinc ion. In addition, theoretical calculations of binding free energy were performed highlighting the key role of specific residues for protein-inhibitor recognition. Overall, these data are very useful for the development of new inhibitors with higher selectivity and efficacy for various CAs.


Introduction
Carbonic anhydrases (CAs) are metalloenzymes that catalyse the reversible hydration of carbon dioxide to bicarbonate and proton 1,2 . CAs are widespread in organisms belonging to all life kingdoms (i.e. bacteria, archaea, and eukarya) and evolved into eight distinct families, namely a, b, c, d, f, g, h, and i. a-CAs are present mainly in vertebrates, fungi, protozoa, corals, algae, in the cytoplasm of green plants, and in some bacteria 3 . b-CAs have been found in bacteria, algae, and chloroplasts of both monocotyledons and dicotyledons, as well as in many fungi and archaea 4 . c-CAs have been reported in archaea, bacteria, and plants 5 , whereas dand f-CAs are found only in marine photosynthetic eukaryotes [6][7][8] . gand h-CAs have been discovered in Plasmodium species and in the marine diatom Phaeodactylum tricornutum, respectively 9,10 . Finally, the new i-CA subclass was recently identified in the marine diatom Thalassiosira pseudonana 11 . All human (h) CAs belong to the a-class, with 12 catalytically active isoforms so far identified. These enzymes are extensively distributed in many tissues and organs, where they are involved in several essential physiological processes such as pH and CO 2 homeostasis, respiration, transport of CO 2 /bicarbonate, electrolyte secretion, biosynthetic reactions, bone resorption, and calcification 2,12 . Consequently, their dysregulated expression level and/or abnormal enzymatic activity may be related to pathological conditions. For this reason, hCAs have been recognised as targets for the design of inhibitors or activators useful for biomedical applications 1,[13][14][15] .
To date, several classes of CA inhibitors (CAIs) have been biochemically and structurally characterised, with sulphonamides and their isosteres (sulphamates and sulphamides) being the most studied 1,16 . However, very often the developed inhibitors were poorly selective, being able to inhibit indiscriminately all or most of the hCA isozymes 1,2,16 . With the aim to obtain more selective CAIs, recent years saw the continuous development and testing of molecules containing new chemotypes such as the carboxylates, polyamines, sulfocoumarins, and coumarins [17][18][19][20][21] . Among them the benzoxaborole derivatives were recently shown to have good inhibitory properties against aand b-CAs from pathogenic fungi and protozoans 22,23 , as well to be effective on different hCA isoforms (Table  1) 24 . In particular, the lead compound benzoxaborole 1 was shown to act as a selective hCA I and hCA II inhibitor and to bind the enzyme active site in its tetrahedral anionic form independently from the pH used in crystallisation experiments ( Figure 1). Interestingly, two different binding modes were observed: in the first one, one of the hydroxyl groups of the benzoxaborole was anchored to the catalytic zinc ion, completing its tetrahedral coordination (Figure 1(B)), whereas, in the second one, the inhibitor was bound to the zinc ion with two of its oxygen atoms, generating a trigonal bipyramidal coordination geometry of zinc (Figure 1(C)) 24 . The substitution pattern at the benzoxaborole ring was deeply investigated too, showing that benzoxaboroles containing urea and thiourea moieties (compounds 2 and 3 in Table 1) were able to discriminate between the different CA isoforms 24 . Based on these interesting data, in this article, we further extend the characterisation of ureido/ thioureido-benzoxaboroles by means of X-ray crystallography and binding free energy calculations.

X-ray crystallography
Crystals of homemade hCA II in complex with 2a, 2b, 3a, and 3b were prepared by the soaking technique 25 . In particular, hCA II native crystals were grown at room temperature using the vapour diffusion hanging drop method. Equal volumes of protein (10 mg/mL in 0.02 M Tris-HCl pH 8.0) and precipitant solution (1.3 M sodium citrate, 0.1 M Tris-HCl, pH 8.5) were mixed and equilibrated against 1 mL reservoir, containing the same precipitant solution. A few enzyme crystals were then transferred in a 2 lL drop of freshly prepared precipitant solution containing 40 mM inhibitor and 10% (v/v) glycerol as cryoprotectant agent. These crystals were kept in the soaking solution overnight for hCA II/3a and hCA II/3b adducts, whereas a very long soaking time was necessary for hCA II/2a and hCA II/2b (1 week for both cases). Crystals were frozen in a gaseous nitrogen stream prior to the diffraction experiment without further manipulations.
A complete X-ray data set was collected at 100 K by a copper rotating anode generator developed by Rigaku and equipped with  Binding of benzoxaborole to hCA II active site (PDB accession code 5JQ0) 24 showing its tetrahedral (B) or trigonal bipyramidal coordination (C) to the catalytic zinc ion.
a Rigaku Saturn CCD detector. Diffraction data were processed and scaled using program HKL2000 (HKL Research, Inc., Charlottesville, VA) 26 . Data collection statistics are reported in Table 2. The initial phases of hCA II/inhibitor structures were calculated using the atomic coordinates of the native enzyme with waters removed (PDB entry 1CA2) 27 . Electron density for all inhibitors was observed in the difference map after a single round of refinement. A model for each inhibitor was then built and introduced into the atomic coordinates for further refinement, which proceeded to convergence with alternating cycles of water addiction, manual rebuilding with the O program 28 , and energy minimisation and B-factor refinement with the CNS program 29,30 . Topology files for all compounds were obtained using the PRODRG server 31 . Refinement statistics for all hCA II/inhibitor adducts are summarised in Table 2. Coordinates and structure factors have been deposited in the Protein Data Bank (accession codes 6RVF, 6RVK, 6RVL, and 6RW1 for hCA II/2a, hCA II/2b, hCA II/3a, and hCA II/3b, respectively).

Computational analysis
Theoretical calculations were carried out on the four crystallographic complexes of hCA II with 2a, 2b, 3a, and 3b as well as on two model complexes with 3a * and 3b * . The latter models were obtained by substituting the oxygen atom of urea moiety with a sulphur atom in the hCA II/2a-b crystal structure, using the Builder module of Insight II package (Insight2000, Accelrys, San Diego, CA). The partial atomic charges for ligands were obtained by quantum mechanical (QM) calculations (B3LYP/6-31G Ã ) using the Gaussian16 software (Gaussian, Inc., Wallingford, CT) 32 through the restrained electrostatic potential (RESP) fitting procedure as implemented in the PyRED server 33,34 . Since benzoxaborole derivatives bind to the catalytic metal ion in their tetrahedral anionic form 24 , the total charge for ligands was set at À1 e. According to our previous works 35 , a charge of þ1.5 e was used for the zinc ion. The General AMBER 36 and the AMBERff14SB force fields 37 were used for the ligands and proteins, respectively 38 . Van der Waals parameters for the Zn 2þ ion were adopted from the work of Li and Merz 39 (r ¼ 1.271; E (kcal/mol) ¼ 0.00330286). Since boron atom is not parametrised in Amber-derived force fields, it was substituted with carbon, as reported in other modelling studies 40,41 , however, retaining the partial atomic charge of boron atom computed by PyRED as described above. The binding free energies (DGbind in kcal/mol) for crystallographic complexes and the theoretical models were calculated using the molecular mechanics/generalised Born surface area (MM/GBSA) method 42,43 implemented in AmberTools18 44 . Moreover, to identify the key protein residues responsible for the ligand-binding process, the binding free energy was decomposed on a per-residue basis.
where Ii(hkl) is the intensity of an observation and < I(hkl)> is the mean value for its unique reflection; summations are over all reflections.
R-free is calculated as for R-work, but from data of the test set that was not used for refinement ( According to MM/GBSA method, the binding free energy was estimated as follows: where DG bind is the binding free energy and G complex , G protein , and G ligand are the free energies of complex, protein, and ligand, respectively. In details: where DG bind is the binding free energy in solution, DE MM is the molecular mechanics energy, which comprises van der Waals (DE vdW ) and electrostatic (DEelec) contributions; DG sol is the solvation energy, and includes electrostatic (DG GB ) and nonpolar (DG SA ) interactions. The electrostatic solvation energy (DG GB ) is evaluated using the Generalised Born method 45 , and the nonpolar contribution is computed by the Linear Combination of Pairwise Overlaps (LCPO) method 46 . TDS is the change of conformational entropy on ligand binding, which is not included in our calculations. Indeed, for comparison of similar ligands, it is acceptable to exclude the entropy contribution in practise 43,47 .

Results
To elucidate the binding mode of ureido/thioureido-benzoxaborole derivatives to the CA active site, we carried out a crystallographic study of the ubiquitous hCA II in complex with compounds 2a, 2b, 3a, and 3b (Table 1), characterised by the presence of a phenyl ring anchored to the urea/thiourea group through a linker of variable length. hCA II was selected as a model isoform for these studies since it is easy to crystallise and numerous data are available on adducts that it forms with inhibitors belonging to different classes 1 . Crystals of hCA II adducts were obtained by soaking experiments following a protocol well described in literature 24,48 . In particular, native crystals of hCA II were grown using the hanging drop vapour diffusion technique at pH 8.5 and the obtained crystals were then transferred into a freshly prepared precipitant solution containing also the inhibitor. Data collection and refinement of all structures were performed as reported in the experimental section.
Analysis of the electron density maps since initial stages of the crystallographic refinement revealed that, as previously reported 24 , all benzoxaborole derivatives bind to the catalytic metal ion in their tetrahedral anionic form. Electron density maps were well defined for the whole inhibitor molecule in the case of hCA II/ thiourea adducts, whereas a greater conformational variability was observed for the phenyl substituent of urea derivatives (Figure 2). The inhibitor binding did not alter hCA II three-dimensional structure. Indeed, the r.m.s.d. values calculated by superposition of all the Ca atoms of the hCA II/inhibitor adducts with those of the native protein were very low (r.m.s.d. values in the range of 0.18-0.21 Å).
A careful analysis of the hCA II/2a and hCA II/2b adduct structures showed that both inhibitors bind to the enzyme active site in a similar way, coordinating the catalytic zinc ion in tetrahedral geometry through one of the hydroxyl groups of the benzoxaborole moiety (Figure 2(A,B)), and being stabilised by many other polar and hydrophobic interactions. The conformations of the two inhibitors in the catalytic cavity are very similar (Figure 3(A)) with few differences related to the lengths of the phenyl-urea tails. The comparison of these structures with the previously reported compound 2c 24 showed that, despite a rather similar orientation of the benzoxaborole skeleton in the active site, the latter compound coordinates the zinc ion in a trigonal bipyramidal geometry (Figure 3(A)). Based on this finding and in agreement with previously reported data 24 , it is reasonable to assume that the two zinc coordination geometries, the tetrahedral and the trigonal bipyramidal, are in this case energetically equivalent and small changes in the orientation of the benzoxaborole skeleton can allow the transition from one to the other.
Superposition of compounds 3a and 3b, when bound to the enzyme active site, showed that, as observed for the urea-containing compounds, the binding modes of the two thiourea derivatives are similar to each other (Figure 3(B)). In this case, the two inhibitors bind the metal ion with two of their oxygen atoms, generating a trigonal bipyramidal coordination geometry. Additional polar and hydrophobic interactions contribute to stabilise the binding (Figures 2(C,D)).
In contrast, great differences are observed when the two groups of compounds are compared with each other. Indeed, even if compounds 2a and 2b differ from 3a and 3b only for one atom (an oxygen atom instead of a sulphur one) (see Table 1), their arrangement in the enzyme active site is completely different (Figure 4). One of the main differences is related to the geometry of the urea/thiourea moiety, which is trans-trans in the case of 2a and 2b and trans-cis in 3a and 3b ( Figure 5). A question comes spontaneously to mind at this point: have urea and thiourea moieties an intrinsic preference for trans-trans and trans-cis geometry, respectively, that leads them to adopt a different binding orientation within the active site or are there specific protein/inhibitor interactions which influence this geometry and are responsible for their different binding mode to the enzyme?
To answer this question, a detailed analysis of urea and thiourea bis-phenyl derivatives contained in the Cambridge Structural Database (CSD) 49 was carried out. This analysis revealed that these molecules do not have intrinsic conformational preferences, but they can assume both trans-trans and trans-cis conformations depending on the nature of the phenyl substituents and on the chemical environment (data not shown). Thus, the different binding geometries assumed by urea and thiourea derivatives within the hCA II active site are induced exclusively by specific interactions with the enzyme residues. Based on the inhibition data reported in Table 1, showing that ureido-benzoxaboroles have greater affinity for hCA II with respect to thiourea derivatives, it can be supposed that the binding conformation observed for compounds 2a and 2b is energetically more favoured with respect to that adopted by compounds 3a and 3b. Thus, to understand why thiourea-containing inhibitors do not assume the same conformation of urea-containing ones, binding free energy calculations were performed by using the MM/GBSA theoretical method 42,43 . This method allows decomposing the binding free energy on a per-residue basis in order to identify key protein residues responsible for the inhibitor binding mode. Calculations were carried out on the four crystallographic complexes as well as on two model complexes (hCA II/3a * and hCA II/3b * ) obtained by substituting the oxygen atom of urea moiety with a sulphur atom in the hCA II/2a and hCA II/2b crystal structures. The latter two structures represent hypothetical models in which thiourea derivatives would adopt the same binding conformations observed for the corresponding urea derivatives. Table 3 reports the protein residues giving a major contribution to the ligand-binding energy in the case of the crystallographic complexes, whereas Table 4 reports the same data for the model complexes. In the latter case, the energy differences with respect to the crystallographic structures are reported. The zinc ion contribution to ligand binding energy is not reported since, as reported in the literature, it is affected by the overestimation of the electrostatic interactions due to the high positive charge of Zn 2þ 50 .
For what concerns the crystallographic complexes, the protein residues mainly contributing to the binding are the same for all the ligands with quite comparable extents. In detail, Phe131, Val135, and Pro202 participate in the inhibitor binding through van der Waals (vdW) interactions (Table 3). These residues form a wide hydrophobic pocket at the mouth of the catalytic cleft, which is able to accommodate the phenyl ring of each ligand (Figure 2). Val143 and Leu198, located in the interior of the binding pocket, are involved in stabilising vdW and electrostatic interactions with the ligands, interacting with the benzoxaborole ring (Table 3). In addition, Thr199 and Thr200, at the bottom of the active site, establish electrostatic interactions with the inhibitors. Finally, Gln92, located at the entry of the active site, also  contributes to inhibitor stabilisation. However, it is worth noting that a very strong polar interaction is observed in the case of hCA II/3b adduct, due to a hydrogen bond between the sulphur atom of thiourea moiety and Gln92 side chain, which is absent in the case of the other derivatives (Table 3 and Figure 2).
Regarding the model complexes (Table 4), the energy differences with respect to the crystal adducts are close to zero for all residues except for Phe131, which shows a destabilising energy contribution equal to 15.0 and 1.5 kcal/mol for derivatives 3a * and 3b * , respectively. These data indicate that thiourea derivatives would experience a significant destabilising interaction with Phe131 if they would adopt the same binding conformation of urea derivatives, likely due to an increase of the steric hindrance between their CS group and the Phe131 ring. These data suggest that the thiourea moiety turns towards a trans-cis conformation, to avoid this unfavourable interaction.
Our energetic calculations can help also to rationalise the interesting inhibition profile of the thiourea derivatives against the different hCA isoforms, with a particular focus on hCA IX and XII which have been recently recognised as valuable targets for cancer treatment and diagnosis 51 . Indeed, inhibition assays showed that these compounds inhibit with greater efficiency the tumourassociated hCA IX and XII with respect to the ubiquitous hCA II (Table 1). This can be explained based on the substitution of Phe131 with a valine and an alanine in hCA IX and XII (Figure 6), respectively. Indeed, these two smaller residues likely do not create the steric hindrance, originated by Phe131 in hCA II, for the positioning of the inhibitor CS group and permit an energetically more favoured arrangement of thiourea derivatives. These hypotheses are in agreement with previously reported structural studies 1,48,52 , which evidenced the important role of CA residue in position 131 in protein-inhibitor recognition and support the idea that region encompassing residues 131-135 constitutes a "hot   Only residues contributing more than À1.0 kcal/mol to the binding are reported. DE vdW : van der Waals contribution; DE elec : electrostatic contribution; DG GB : generalised-Born solvation contribution; DG SA : non-polar solvation contribution. zone" 53 to be targeted for the design of novel CAIs with increased selectivity.
Altogether, studies reported here demonstrate that ureido-and thioureido-benzoxaborole derivatives are an interesting class of versatile CAIs. Indeed, these molecules can bind the enzyme active site assuming different binding conformations and coordination geometries depending on the interactions established with the enzyme active site residues. This feature can be efficiently used to design molecules with enhanced selectivity for pharmaceutically relevant hCA isoforms.